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We explore the long-time dynamics of a system of identical charged particles trapped on a closed 
helix. This system has recently been found to exhibit an unconventional deformation of the linear 
spectrum when tuning the helix radius. Here we show that the same geometrical parameter can 
affect significantly also the dynamical behaviour of an initially broad excitation for long times. In 
particular, for small values of the radius, the excitation disperses into the whole crystal whereas 
within a specific narrow regime of larger radii the excitation self-focuses, assuming finally a localized 
form. Beyond this regime, the excitation defocuses and the dispersion gradually increases again. 
We analyze this geometrically controlled nonlinear behaviour using an effective discrete nonlinear 
Schrodinger model, which allows us among others to identify a number of breather-like excitations. 


PACS numbers: 37.10.Ty, 37.90.+j, 45.90.+t, 05.45.-a 

I. INTRODUCTION 

Whereas the harmonic approximation of interactions 
provides valuable information about the stability and the 
propagation of small amplitude excitations in crystals 
formed by interacting particles, their real-time dynamics 
as well as their thermal and transport properties are typ¬ 
ically subject to some degree of nonlinearity pQ. Among 
the most prominent manifestations of such a nonlinearity, 
are the self-focusing or self-trapping mu of initial wave 
packet excitations and the existence of non-spreading ex¬ 
citations like breathers and kinks For discrete sys¬ 

tems, a prototype equation incorporating these features 
is the so-called discrete nonlinear Schrodinger (DNLS) 
equation consisting of a linear (dispersive) coupling and 
a cubic nonlinear term |8j, used to model plenty of sys¬ 
tems ranging from coupled optical waveguides mm and 
Bose-Einstein condensates Q2HII] to transport in DNA 
molecules [HHUl- The standard spatial arrangement of 
sites in most of such one dimensional (ID) studies is that 
of a straight equidistant chain in which the coupling (hop¬ 
ping) is restricted to nearest-neighbors (NN). 

Non-trivial lattice geometries for ID discrete nonlin¬ 
ear systems have also been studied and have been found 
to lead to intriguing new phenomena owing to the in¬ 
terplay between geometry and nonlinearity. In particu¬ 
lar, in curved ID lattices embedded in a 2D space the 
bending can act as a trap of excitations and induce a 
symmetry breaking of nontopological solitons [18]. For 
a 3D space, the helicoidal lattice structure, such as that 
of DNA molecules, is found to enhance the existence and 
stability of discrete breathers m 20 . Furthermore, a 
curved geometry has been proven to induce nonlinear¬ 
ity in systems where the underlying interactions are har¬ 
monic [2T1 [22] . 

In the present work we examine the interplay be¬ 
tween nonlinearity and geometry in a system of identical 
charged particles, confined on a curved ID manifold em¬ 


bedded in the 3D space, namely a closed (toroidal) helix. 
In a previous work [23], we have shown that in such a 
system, a tuning of the geometry controlled by the he¬ 
lix radius, leads to an unconventional deformation of the 
phononic band structure including a regime of strong de¬ 
generacy. As a consequence, the propagation of small 
amplitude localized excitations is affected significantly 
and a specific geometry exists at which the excitations 
remain localized up to long times. A natural question 
therefore arises, what would be the long time dynamics 
of a general excitation and whether there is some geomet¬ 
rically controllable degree of nonlinearity inherent in the 
system which can alter the propagation characteristics. 

We provide an answer to this question by studying the 
time propagation of an initially broad excitation on the 
crystal of charges. We find that for values of the he¬ 
lix radius far from the degeneracy regime, the excitation 
initially spreads with multiple subsequent revivals due 
to the closed shape of the crystal. Within the degener¬ 
acy regime, however, the initial excitation focuses in the 
course of propagation reaching finally a rather localized 
state, serving as a hallmark of the existing nonlinearity. 
In order to quantify this nonlinearity we construct an ef¬ 
fective DNLS model with additional nonlocal nonlinear 
terms [24] [25]. Such a model is found to capture qual¬ 
itatively well the localization and dispersion features of 
the original dynamics, providing a deeper insight into the 
observed effect. Even more, it gives us the opportunity 
to identify some discrete breather-like excitations at the 
degenerate geometry, thus adding to the dynamical pic¬ 
ture. 

The structure of this work is as follows. In Sec. II 
we describe our system commenting also on its linearized 
behaviour. In Sec. Ill we present our results for the 
time evolution of an excitation in the crystal for different 
geometries. In Sec. IV we construct an effective DNLS 
model for our system at the geometries of interest and 
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in Sec. V we use it to identify some breather-like excita¬ 
tions. Finally Sec. VI contains our conclusions. 

II. SETUP AND LINEARIZATION 

We consider a system of N identical charges of mass 
m o which interact via repulsive Coulomb interactions and 
are confined to move on a ID toroidal helix, parametrized 
as 


( (R + r cos(u)) cos(mx)\ 

(R-h r cos(u)) sin(au) J , uE [0, 2M7r]. (1) 
rsin(u) J 

In eq. 0 R denotes the major radius of the torus 
( Fi g0a)), h is the helix pitch and r the radius of the 
helix (minor radius of the torus), whereas a = stands 
for the inverse number of windings M = The to¬ 

tal effective interaction potential which results from the 
constrained motion of the charges on the helical man- 
ifold reads V{u u u 2 ,... u N ) = \ E^=i Mi |r(^)-rK)| > 
where Uj denotes the coordinate of particle j and A is the 
coupling constant characterizing the standard Coulomb 
interactions. 

We choose dimensionless units by scaling all our physi¬ 
cal quantities (e.g. position x, time t and energy E) with 
A, m o and 2h/'K as follows 


xi r ~ / \i r 3 ~ 2 Eh „ ~ 

*=2S’ ' = t Vw?' m ° = ’ x= • 


omitting in the following the tilde for simplicity. 

At commensurate fillings, i.e. M = nN, n = 1,2,... 
with the filling factor being z^ = l/n<l,itis found that 
for values of the helix radius r up to a critical point r c 
the ground state configuration of such a system is the 
equidistant polygonic configuration = 2 (j — 1) irn 
(Fig. 0 a)). This configuration loses its stability at r c 
undergoing a zig-zag bifurcation [23]. 

We focus in this work on the dynamical behaviour of 
charged particles confined on the toroidal helix in the 
region r < r c , where the ground state is still the poly¬ 
gonic one. We have shown in [23] that in such a region 
the linear spectrum of the system changes dramatically 
with tuning the radius of the helix r, a fact that crucially 
affects the propagation of small amplitude localized ex¬ 
citations. Specifically, it was found that the width of the 
linear spectrum decreases as one approaches a point r<i of 
strong degeneracy from below and increases again beyond 
that point, while interchanging the character between the 
eigenmodes corresponding to the highest and the lowest 
frequencies (Fig. [ljb)). In fact, since the degeneracy is 
not complete, it is better to refer to a degeneracy regime 
within which the inversion of the spectrum is gradually 
achieved while its width remains small (Fig. 0b) (in¬ 
set)). We consider in this work six different geometries, 



FIG. 1. (color online), (a) Equidistant configuration of ions 
confined on the toroidal helix for v — \ and N — 8. The yel¬ 
low arrows indicate the initial velocities of the particles, (b) 
Highest (solid blue line) and lowest (solid red line) frequencies 
of the linearization spectrum around the equilibrium equidis¬ 
tant configuration as a function of r for N = 60 particles. 
The black dots and the empty diamonds refer to the frequen¬ 
cies corresponding to the center of mass and the out of phase 
mode respectively. The vertical lines mark the radii of the 
helix we use in our calculations. The small insets depict the 
form of the respective vibrational band structures uo(k) at the 
corresponding values of r. (c) Initial local energy E n profile 
as a function of the particle index n. 


each corresponding to a different value of r, covering all 
the regions with a qualitatively different linear spectrum 
(Fig. [ljb)) from the ring limit (ro = 0) to the degeneracy 
(^ 3 ,^ 4 ) and the inversion (r 5 ) regime. We focus on the 
case of half-filling v — ^ for N = 60 particles. 

III. TIME PROPAGATION OF A GAUSSIAN 
EXCITATION 

In this section we present and discuss the dynamical 
response of our system to an initial excitation. Although 
the physical results are in principle independent of the 
exact character of this excitation and the means used for 
its quantification, the determination of both is essential 
for the illustration and the theoretical description of our 
findings. 

Dealing with classical systems and seeking an excita¬ 
tion measure whose total amount is conserved in time, 
the natural choice is a, to be defined, energy distribution 
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associated with each particle, referred to hereafter as lo¬ 
cal energy E n . Whereas the kinetic energy K consists of 
parts allocated to each individual particle, the potential 
energy cannot be uniquely partitioned, yielding different 
definitions of local energy [261428] . Aiming for them to 
be strictly positive for all possible excitations (a consid¬ 
erably non-trivial requirement for systems with Coulomb 
interactions), we define our local energies E n in a rather 
unconventional way, focusing on a positive decomposi¬ 
tion of the harmonic interaction term [28] • Our complete 
definition and a more detailed discussion of local energies 
are provided in the Appendix. 

We start with an initially broad excitation of a Gaus¬ 
sian profile in terms of local energies (Fig. [ljc)) . Since 
the local energies E n depend trivially on the particles’ 
velocities (contrary to what is the case for the particles’ 
positions), the most straightforward way to obtain such 
a Gaussian local energy profile is by exciting the parti¬ 
cles with a suitable velocity distribution. Of course, for 
a given local energy profile the magnitude of such a ve¬ 
locity distribution can be uniquely determined, but there 
is a freedom in the direction of the velocity for each par¬ 
ticle. We choose here all the velocities to point in the 
positive direction (Fig. B a ))- 

The ensuing dynamics is shown in (Fig. [2](a)-(f )). Ob¬ 
viously, the time evolution of the Gaussian excitation 
possesses a drastic dependence on the geometry, con¬ 
trolled by the helix radius r. In particular, for r < 
(Fig. [2](a) ,(b)) the excitation spreads into the whole crys¬ 
tal and refocuses almost periodically at the time instants 
when the left and the right propagating parts of the ex¬ 
citation meet and superimpose at the diametrically op¬ 
posite point of the closed helix. As discussed in [23] the 
spreading velocity decreases as r is increased, following 
the width of the linear spectrum (Fig. 0b)). As the 
width becomes smaller the features of the time evolution 
alter significantly. Already at the point the excitation 
does not spread any more into the crystal, but it alter¬ 
nately focuses and again defocuses to its initial shape 
(Fig- 0 c)). Even more surprisingly, within the degener¬ 
acy region (Fig. 0d) ,(e)) the initial excitation under¬ 
goes a focusing after some time scale tp (tp « 4000 
for (d), tp ~ 6000 for (e)). Subsequently, the wave 
packet loses its smooth envelope and fragments into a 
number of highly localized excitations. Depending on r, 
the routes towards such a localized state can be different, 
with the wave packet evolving initially one central peak 
(Fig. [2](d) ) or two side peaks (Fig. [2|e)). Another inter¬ 
esting feature within the degeneracy regime is that the 
reflection symmetry of the initial excitation profile can 
break in the course of propagation, attaining after some 
time a significantly asymmetric form (Fig. [2](d) ). The di¬ 
rection of the asymmetry depends on the direction of the 
initial velocities of the particles. Beyond the degeneracy 
(Fig. [2|f)) the spreading of the excitation into the crystal 
reappears with a periodic refocusing but the propagation 


pattern is much different owing to the inverted form of 
the vibrational band structure. 

In order to quantify the degree of focusing or localiza¬ 
tion of the excitation, we examine the time dependence 
of the normalized participation ratio 


V = 


1 (£n= lE n ) 2 
N £n=l E n 


( 2 ) 


Evidently, this quantity can take values between 1/N 
and 1, with V = 1 signifying the case of a completely 
extended excitation where the energy is equipartitioned 
between all the particles and V = 1/N marking the op¬ 
posite case of a fully localized excitation in a single par¬ 
ticle. Note that the local energies E n in (J2| should be 
non-negative for the definition to make sense. 

Our results, presented in Figs. |2jg)-(l), support our 
discussion above. Especially the focusing of the excita¬ 
tion after tp is evident in Figs. [2jj)-(k). However, the 
subsequent drop in V is much stronger in Fig. 0j) than 
in Fig. 0k), in line with the observation that at r = r% 
the final localized state consists of less excited particles 
(Fig. |2jd)) than at r = (Fig. (2je)). 

In [23] it was demonstrated that a small localized ini¬ 
tial excitation does not spread significantly for short 
times at the degeneracy point, in contrast to the be¬ 
haviour for other geometries. This fact was understood 
solely by an inspection of the linearization spectrum. 
Here, however, the situation is different. Not only the 
complete absence of spreading, but especially the exis¬ 
tence of self-focusing calls for an account of the under¬ 
lying nonlinearity. This is further emphasized and sup¬ 
ported by Figs. |2jg)-(l) where the results of the propaga¬ 
tion within the harmonic approximation of the potential 
are also displayed. As long as the total amplitude of the 
initial excitation is small enough, the harmonic approx¬ 
imation works well. As the amplitude is increased this 
approximation will start to fail, and nonlinear effects are 
expected to show up. The results of the present work 
suggest that except for the amplitude, also the geome¬ 
try allows to control the importance of the nonlinearity. 
Specifically, for the given amplitude and for geometries 
far from the degeneracy regime (Figs. 0g)>(h),(l)) the 
harmonic approximation qualitatively reproduces the ex¬ 
act time evolution of the participation ratio, although, 
as should be expected, there are quantitative deviations. 
In contrast, close to and within the degeneracy regime 
the harmonic approximation fails completely, predicting 
a spreading and an extended form of the excitation, in¬ 
stead of localization. This makes it clear that regard¬ 
ing the focusing, we indeed encounter a nonlinear phe¬ 
nomenon. 

Complementary information about the spreading or 
localization of the Gaussian excitation can be obtained 
from its time evolution in the reciprocal space, i.e. in 
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FIG. 2. (color online) (a)-(f) Time evolution of the initial Gaussian excitation presented in Fig. 0 c ) for N = 60, v = 1/2 
for increasing r corresponding to the points (a) r*o, (b) t*i, (c) 7 * 2 , (d) 7 * 3 , (e) 7*4 and (f) 7*5 marked in Figjljb). Colors encode 
the values of local energy E n for each particle n and time t. For the same values of r panels (g)-(l) depict the time evolution 
of the normalized participation ratio V of the excitation. The solid blue lines are the results from our numerical simulations 
corresponding to (a)-(f), whereas the dashed red lines correspond to the results for a harmonic approximation of the potential. 


terms of the wave numbers k G [—^, , where As de¬ 

notes the arc length interparticle distance of the equidis¬ 
tant ground state configuration. To this extent, we ex¬ 
amine how the discrete Fourier transform of the local 
energy profile evolves in the reciprocal space (Figs. [3^a)- 
(f)). The initial excitation, being overall extended in the 
coordinate space (Fig. [ljc)), appears to be rather local¬ 
ized around the wave vector k = 0 in the k- space. It 
remains localized as expected for r < (Figs. !a),(b)), 
but when approaching r& more wave vectors in the vicin¬ 
ity of k = 0 become excited (Fig. |3|c)). At the degen¬ 
erate geometries (Figs. Id) ,(e)) the excitation expands 
rapidly in the reciprocal space, populating after the char¬ 
acteristic time tp almost all the wave numbers k. The 


space localized solution consists therefore of most of the 
k modes, with the population of the initial k = 0 mode 
being dominant. At r > (Figs. |3jf)) the excitation re¬ 
mains, as for small r, in the narrow vicinity of the k = 0 
mode. 

Before proceeding with our study of the nonlinear be¬ 
haviour, let us note that for the results presented here we 
have used a rather small initial excitation with a total en¬ 
ergy of the order of 1% of the ground state energy per 
particle Eqs/N. For larger amplitudes the self-focusing 
can occur also for smaller radii, i.e. at a greater distance 
from the degeneracy region. 

IV. EFFECTIVE NONLINEAR MODEL 

The dynamics analyzed in the previous section is char- 
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FIG. 3. (color online) Absolute value of the discrete Fourier 
transform of the local energy excitation profile as a function 
of time and wave number k. Subfigures (a) to (f) correspond 
to increasing r, ranging from ro to 7*5 as indicated in Fig. [l] 

(b). 


acterized by a self-focusing process of excitations for the 
case of degenerate geometries and therefore suggests a 
prominent role of the nonlinearity. We aim in this sec¬ 
tion to identify and quantify the leading nonlinear terms 
as well as to derive a DNLS effective model in the region 
of degeneracy. As a result we will, among others, gain 
insight into the excitation amplitudes and the time scale 
tp for localization. 

A. The dominant nonlinear terms 


Since the initial excitation is small enough we attempt 
to identify the dominant anharmonic terms by expanding 
the potential around the equilibrium configuration } 
up to fourth order. It is advantageous to do so in the arc 
length parametrization s so that the final Euler-Lagrange 
equations of motion and particularly the kinetic terms 
would assume the standard form [29] . To this purpose 
we calculate the matrices 


rr _ &V 

u dsidsj 




, Gij - 


d 3 V 


lJ dsjdsj 




(3) 


corresponding to the Hessian and the matrix of the third 
derivatives respectively. Terms involving more than two 
different positions are equal to zero, since the total poten¬ 
tial V is a sum of exclusively two-body potential terms. 

For the fourth derivative terms we define two further 
matrices as 


d A V 

ij ~ ds?ds* 


5 Qij 


d A V 




ds^dsn 




(4) 


for i 7 ^ j and 


separating the derivatives of the same order in Si, Sj from 
those with a different order and splitting the diagonal 
derivative terms (involving differentiation with respect 
to a single position) in two. The calculation of the above 
derivatives in the arc length parametrization can be car¬ 
ried out by using the respective derivatives in the u coor¬ 
dinate space as well as the known relations for derivatives 
of inverse functions. Under these considerations and de¬ 
noting Sj — = Xj the potential reads 


^ 1 d 4 V 

M u -Qu- 2 


V « E gs + XiHijXj + 1 x?GijXj 

i,j =1 i,j =1 

1 JV J AT 

+ 24 E x i M a x2 j + ^E x iQij x j ( 5 ) 

i,j=l i,j =1 

leading to the equations of motion 

X n — ^ ^ H n jXj H nn X n — ^ ^ X n G n jXj 


j^n 


j^n 


— x 2 --\^ x m x 2 -- M x 3 

~ £ / J •^n lv - L nj^j ^ lv - L nn a ^m 


j^n 


j^n 


j^n 

£ X nQnjXj — — Q 


6 jn x j 5 


(6) 


' j/n 


which consist of harmonic, quadratic and third order non¬ 
linear terms. All the matrices appearing in these equa¬ 
tions are symmetric except for the matrix Gij, relating 
to the quadratic force terms, which is fully antisymmet¬ 
ric. This can affect significantly the symmetry of the ex¬ 
pected solutions. If the quadratic nonlinear force terms 
are ignored then the equations of motion § possess the 
symmetry x n —>• — x n which also results in symmetric ex¬ 
citations keeping their symmetry in the course of prop¬ 
agation. The quadratic force terms, however, break the 
reflection symmetry and allow for an asymmetric evolu¬ 
tion of initially symmetric excitations as the one observed 
in Fig. §d). 

Let us mention at this point that in the ring limit r = 
ro = 0 , where a separation of the center of mass holds, 
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the matrices involved in eq. © are not independent but 
obey the relations 


_^ ^ T 3 

H nn ^ ^ H n j , M nn — ^ ^ Qnji Mnj 

j^n j^n 


yielding 


X n — H n j (xj X n ) Gjn(Xj ■ 




j/n 


1 

24 


Qnj{ x j x n) (J) 

j^n 


which for only nearest neighbor (NN) couplings leads to a 
Fermi-Pasta-Ulam kind of equations of motion WL with 
both quadratic and cubic nonlinear interactions. 

Before proceeding, we note that the matrix elements 
depend, due to the symmetry of the 
ground state configuration, only on the index difference 
m = i — j. Therefore, when referring to these elements in 
the following we will use the notation H m , G m , M m , Q m . 


B. The DNLS model 


The effect of localization of the initial wave packet 
occurs, as we have observed in Sec. Ill, in the regime 
close to degeneracy. There, the diagonal terms of the 
Hessian provide the dominant contribution to the linear 
spectrum as the off-diagonal ones are very small. Since 
\Hq\ >> 2|Hi| (Fig. J^a)), one can use in that regime 
the so-called rotating wave approximation (RWA) [51] . 
assuming that the position coordinate can be described 
as 



FIG. 4. (color online), (a) The ratio between the diagonal 
and off diagonal elements of the Hessian 2 \hI\ as a function 
of r. (b)-(e) Coefficients Ai, B,Ci, Di defined in eq. (10) as a 
function of r. The solid blue lines depict the coefficient values 
for 1 = 1, whereas the black doted lines for l = 2 and the red 
dashed ones for l = 3. The cyan (shaded) area marks the 
region of applicability of the RWA as derived from the panel 
(a). The vertical lines correspond to the values of r we study 
with the DNLS model. Note that three of them namely r 2 ,r *3 
and r 4 are the same as those studied in the previous section 
(Fig. [ljb)) whereas 7*45 lies between 7*4 and 7 * 5 . The insets 
provide a zoom into the region of interest. 


© acquire the form 


Xn(t) = + K(t)e iMot , (8) 

with T n (£) a slowly varying amplitude and cUq = Ho 
yielding a fast oscillating phase e ±luJot . Apart from 
the requirement of a weak dispersion, a condition for 
a sufficiently weak nonlinearity has also to be satisfied 
[511 | , namely \Hq\ | |M 0 | max(a^(0)) 2 with max(a^(0)) 
the maximum initial displacement (or respectively mo¬ 
mentum) of a single particle. This criterion is satis¬ 
fied as well in our case, since the initial conditions we 
have used lead to \Hq\ of about 100 times larger than 
| |Mo| max(a^(0)) 2 . 

Using the ansatz (|8| and making the assumption 
< cjo 14/ n |, as well as neglecting the rapidly oscil¬ 
lating terms with frequency 3u;o5 the equations of motion 


N/2 


= Y M^n-l + *n+j) + B |^„| 2 ^ n 

) + + *n+l) 


2«' n (i«' n _ i | 2 + i^ n+ y 


1=1 

N/2 

1=1 
N/2 

+ Y D ' [ 2 l^l 2 (*»-« + *«+l) + (K-l + K+l) 

1 = 1 

+ + \^ n+l \ 2 ^ n+l 


N/2 


Y Fi [(2tt„ (*„+, - _ tf n+l + * 2 _ ; 


o -iuj 0 t 


1=1 


+ (2tf n - K-i) ~ I*n+«| 2 + l^n-;| 2 ) e iwot 


(9) 
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with 


A = — R — Ml r = Ml 

1 2w 0 ’ 2 w 0 ’ 1 12w 0 

n = ^ zr = 

1 16w 0 ’ ‘ 12w 0 ' 


( 10 ) 


The last term in eq. ([9]) comes from the quadratic non- 
linear force term which is responsible for generating sec¬ 
ond harmonics and breaking the reflection symmetry of 
the equations. Within the RWA this term is also ne¬ 
glected due to its fast oscillation in time in comparison 
to the other slowly varying terms. Of course such an 
approximation limits our consideration to the general fo- 
cusing/defocusing behaviour of the wave packet excita¬ 
tion, not accounting for secondary propagation features 
such as the exact shape of the excitation pulse during the 
time evolution. 

The actual numerical values of the remaining coeffi¬ 
cients Ai,B,Ci,Di for our system are shown in Figs. 
0b)-( e) as a function of r. All the coefficients in the 
region of validity of our approximation have a definite 
sign except for Ai which changes sign within the degen¬ 
eracy regime. It is also obvious that an increasing index 
l leads to coefficients closer and closer to zero, with even 
the coefficients of / = 2 being already very small. Since 
the major contribution stems from l = 1, we can proceed 
performing a NN approximation which yields 


= i4(¥„_1 + $„+!) + B |^„| 2 T r 


+ c 

+ D 


2^, 


(l^.il 2 + |* n+1 | 2 ) + + *n+l) 


2 |M/„| 2 ($„_1 + *n+l) + K {K-l + *JL+l) 


+ l^n-ll 2 T n _! + |^„+l| 2 T n+ i 


( 11 ) 


where we have omitted the index 1 from the coefficients 
for clarity. Equation can be recognized as a DNLS 
equation with additional nonlinear couplings studied in 
the literature [24, [25] and will be the point of reference 
in the following discussion. 

In order to check the validity of our approximation we 
propagate the Gaussian excitation illustrated in Fig. [ljc) 
according to eq. (11). Note that the quantity T is by def¬ 
inition complex, with its real part at t = 0 relating to the 
displacement Re(T n (0)) = eA91 ( e q. (jg|) an q its imagi¬ 
nary part to the momentum Im(T n (0)) = (eq. ^ 

within the RWA), both in the arc length parametriza- 
tion. Furthermore, the local energy E n = + p^) 

can be expressed as E n = 2cjq |4/ n | 2 , the conservation 
of the total energy being thus linked directly with the 
conservation of A comparison between Figs. 

|5fa)-(c) and Figs. |2][c)-(e) makes it clear that the NN 
DNLS model captures qualitatively very well the features 
of the exact excitation propagation, exhibiting focusing 
and defocusing at r 2 (Fig. |5ja)) and localization after 


a certain time for r within the degeneracy region (Figs. 
|5jb),(c)). Beyond the point 7 * 4 , but still very close to the 
degeneracy (Fig. Qa)), the excitation, although keeping 
its shape for a long time, eventually disperses into the 
whole crystal (Fig. id)). Note that the profile of the 
excitation remains reflection symmetric (with respect to 
the central particle) in the course of propagation for all 
geometries (Figs. |5ja)-(d)), in contrast to what is ob¬ 
served in (Fig. id)), which we attribute primarily to 
the neglected quadratic nonlinear term. 



FIG. 5 . (color online), (a)-(d) Time propagation of the ini¬ 
tial Gaussian excitation of Fig. [TJc) under the DNLS equation 
( 11 ) for increasing r corresponding to the points (a) 7*2, (b) 
7*3, (c) 7*4 and (d) r 45 marked in Fig.j 4 Ta). Colors encode the 
values of local energy E n — 2 uq iT^Mor each particle n and 
time t. 


Most of these results can already be inferred from an 
inspection of eq. ( pTj ) and the parameter A(r) as shown 
in Fig. ^b). When increasing the radius, A changes its 
sign from negative (same sign as B ) to positive. Cor¬ 
respondingly, with increasing r the effective nonlinearity 
changes from attractive (hopping A and nonlinearity B 
of the same sign), leading to self-focusing dynamics (Figs. 
|5ja)-(c)), to repulsive (hopping A and nonlinearity B of 
opposite sign), leading to defocusing (Fig. j5jd)). 

The change of the magnitude of the parameter A af¬ 
fects also the degree and the time-scale of the wave- 






















0 


(xlO- 3 ) 


3.5 


packet localization. A small hopping term is expected to 
slow down the dynamical behaviour of the system since 
it accounts for reduced mobility. A decrease in its magni¬ 
tude leads also to an enhancement of the general nonlin¬ 
ear behaviour since it shifts the weight of the dynamics 
to the nonlinear terms. On a more formal level, these 
statements can be justified within a multiple scale analy¬ 
sis [32]. In fact, one can arrive at the RWA by assuming 
that all the terms of eq. are much smaller than the 
Hnn x n term. Introducing the small parameter e = |A|, 
due to the weak coupling and weak nonlinearity condi¬ 
tions a slow time scale can be defined t = et = | A\ £, which 
is characteristic for the time evolution of T n . A further 

rescaling of T n with y ]^| leads to the standard form 


of our DNLS 

c 


with A = sgn(A), B = sgn(^), 

C = D = and 4/ n = h 1 which the strength 

of the nonlinearity is primarily controlled by the ampli¬ 
tude of the initial excitation T n (0). In the degeneracy 

regime, due to the factor y/|^|, this is about one order 
of magnitude larger than its value for other parameter 
values of r resulting in the observed strong localization 
of the excitation. In the same regime the time scale of 
the time evolution becomes very large, scaling as ^ and 


yielding in our case a characteristic time tp of the order 
of 10 3 . 

Before closing this section, let us note that in the limit 
\A\ —>> 0 the NN neighbor approximation gradually fails 
(Fig. Qb)) and the dispersion coefficients Ai from all the 
neighbors have to be taken into account. Furthermore, 
at long times, after the localization the weak nonlinearity 
condition starts to fail (since the maximum amplitude 
of the excitation increases substantially) making more 
prominent the existence of higher harmonics and leading 
to effects such as an asymmetric time propagation. 


V. BREATHER-LIKE EXCITATIONS 


Having constructed an effective DNLS model we can 
finally pose the question of the existence of breather-like 
solutions in our system in the regime of degeneracy. If the 
coefficients C, D of eq. were zero, our system could 
be described within the so-called anticontinuum limit of 
weakly coupled oscillators, since the linear coupling is al¬ 
ready very small. In this limit, it has been proven that 
discrete breather solutions exist if there is a substantial 
degree of anharmonicity and no resonances with the lin¬ 
ear spectrum [33] . In our case with C,D / 0 one can 
still search for breathers provided, as well, that their fre¬ 
quency does not belong to the phonon spectrum [24]. To 
this end the frequency in the ansatz of eq. (| 8 | should be 
changed from uoo to uj < ujq out of the phonon band, i.e. 

x n (t) = ^nit)e- iM + K(t)e iut . ( 12 ) 


Repeating the steps of the previous section assuming 
that still the RWA holds, we arrive at eq. (13) with an 
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FIG. 6 . (color online). (a)-(d) Time propagation of breather¬ 
like excitations in the degeneracy regime, at r — 7*3 in (a), (b), 
(d) and r — r 4 in (c). The color plots of the first column depict 
the values of the local energy E n whereas the second column 
provides the initial profiles of the excitation in position co¬ 
ordinates, i.e. Au n — u n — Un ^ as a function of the particle 
index n. The third column illustrates the time evolution of 
the participation ratio V for the degenerate geometries and 
the initial conditions of (a)-(d). In Figs, (a)-(d) the initial 
excitation comes only from the displacement A u n of particles 
from their equilibrium. 


additional local term A4A 


= 

+ C' 
+ D' 


A \^ n -1 + tf n+ i) + A^ n + B' \^ n \ 2 


2’f, 


(|^ n _l | 2 + |^„+l| 2 ) + + ^ 2 + i) 


2 |* n | 2 (*„_! + v n+l ) + * 2 n (*u + <p; +1 ) 


+ |^„_i| 2 1 + |^„+i| 2 


(13) 


where A = (w§ - w 2 ) /2w, A! = A^, B' = Bf, 
C' = C— and D f = D-. In the anticontinuum limit 

uj uj 

A' , C", D' = 0 this equation can be easily checked to have 
the solutions T n E {0, y^e z ^ n } with arbitrary phases 
ct> n E M. Simple real breather solutions can thus be 
constructed iteratively as strings of the elements 0 and 

±yStarting from such a breather solution and fol¬ 
lowing the solution trajectories in the parameter space 
A\C\D' with a Newton method |24j (34] we can find 
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possible breather solutions for different values of r corre¬ 
sponding to the geometries of interest. The natural step 
then is to transform back to the x n ,p n coordinates and 
use these as an initial excitation in our system, simulat¬ 
ing their time evolution within the full model equations 
as done in section III. 

Our results are presented in Fig. [6] for different ini¬ 
tial conditions at the degenerate geometries. Evidently, 
many of them (Figs. [6^a)-(c)) keep their solitonic char¬ 
acter in the presence of the full Coulomb interactions, 
a fact that can be justified also by inspecting the time 
evolution of the participation ratio. Indeed, the local 
energy profile of the excitations changes negligibly with 
oscillating potential and kinetic energy parts, mapping to 
oscillating displacements A u n and momenta p n . These 
breather-like excitations include a single particle excita¬ 
tion (Fig. [6ja)) as well as two-particle excitations with 
both opposite (Fig. [6^b)) and equal displacements (Fig. 
§<•))• 

Surprisingly enough, we have observed that the excited 
particles should always be separated by a distance (Figs. 
§b) ,(c)) for the solitonic character to persist, otherwise, 
if neighboring particles are excited there is always a small 
dispersion and a part of the excitation that focuses on 
a single particle, acquiring finally a rather asymmetric 
profile (Fig. [6^d)). This picture is also valid within the 
DNLS framework where such solutions are found to un¬ 
dergo a spontaneous symmetry breaking in the course of 
propagation. 

VI. SUMMARY AND CONCLUSIONS 

We have shown that a system of charged particles con¬ 
fined on a toroidal helix can react in qualitatively differ¬ 
ent ways when exposed to an initial excitation, depending 
on the geometric properties of the confining manifold. In 
particular, while dispersion is the major feature of the 
dynamics for very small helix radii r, the excitation self¬ 
focuses and localizes for values of r in the so-called de¬ 
generacy regime, where the linear coupling is very small. 
Beyond this regime, the time evolution of the excitation 
is characterized by a defocusing which is gradually dom¬ 
inated again by dispersion. Especially the self-focusing 
of the excitation observed for the degenerate geometries 
constitutes a hallmark of the existing nonlinearity in the 
system. Interestingly enough, the nonlinear part of the 
interactions (both local and nonlocal) does not approach 
zero in this regime, nor does it undergo a change in sign 
with increasing r, contrary to the behaviour of the re¬ 
spective linear part. It is primarily this contrast between 
the robustness of the nonlinearity and the variability of 
the linear part which induces the variety of dynamics in 
the system, allowing for its control through the tuning of 
the underlying geometry. 

Furthermore, we have identified the character of the 
leading nonlinear terms and constructed an effective dis¬ 


crete nonlinear Schrodinger model with additional non¬ 
linear couplings which has allowed us to predict and in¬ 
terpret the different responses of the helical chain to its 
excitation. Through this model it has also been possi¬ 
ble to identify some special breather-like excitations in 
our system which propagate in time keeping their shape. 
Overall, we emphasize that the present setup offers many 
possibilities for the implementation of various nonlinear 
models ranging from the Fermi-Pasta-Ulam to different 
discrete nonlinear Schrodinger models with both attrac¬ 
tive and repulsive nonlinearities only by the tuning of 
a single geometrical parameter, i.e. the radius r of the 
helix. 

Regarding its possible experimental realization, it re¬ 
lies on the challenging task of constructing a helical 
trap for charged particles. Although such a construc¬ 
tion might not be straightforward, there have been cer¬ 
tain advances towards this direction in different contexts, 
such as the realization of free-standing helical nanostruc¬ 
tures [35H37] and the design of helical traps for neutral 
atoms [313 [39]. Moreover, our results could be of rele¬ 
vance in studies of energy and excitation transfer in he¬ 
lical molecules, such as DNA. Finally, macroscopic real¬ 
izations of our system could also be possible by using for 
example charged beads as done for the study of polymers 

mm- 

APPENDIX: DEFINITION OF LOCAL 
ENERGIES 

As briefly mentioned in Sec. Ill the definition of lo¬ 
cal energies E n in a system of interacting^particles is not 
unique. For a general potential V = \ Yi ■ ^ . W(ui,uj) 
the most common partition consists in equally distribut¬ 
ing each part of the potential energy among the particles 
involved 

1 N 

E n = x F W ^ U ") + K n ~ Ei n\ (14) 

i,i^n 

with E = \ Y^i^n W ( u i°^ > ) the local ener gy of the 
equilibrium configuration and K n the kinetic en¬ 

ergy of the particle n. Although such a definition is very 
useful in describing systems with nearest neighbor inter¬ 
actions such as oscillator chains EH HU, it can lead to 
problems concerning the description of long-range inter¬ 
acting systems and particularly those involving Coulomb 
interactions, as the one we consider in this paper. The 
major problem is that it generally leads to both positive 
and negative values of E n and thus demands a handling 
of both positive and negative excitations. Negative lo¬ 
cal energies should be avoided if we wish to interpret the 
local energy as the amount or probability of the corre¬ 
sponding particle being excited and proceed in defining 
quantities such as the participation ratio V of eq. 
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FIG. 7. (color online), (a) A schematic illustration of three 
charges on a ring at equilibrium (1,2,3) and after particle 
1 has been displaced (l 7 ,2,3). The dashed lines denote the 
Euclidean distances in equilibrium, whereas the solid ones, the 
distances after the displacement. (b),(c) The time evolution of 
a Gaussian excitation for r = rq. (d),(e) The time evolution 
of a Gaussian excitation for r = rq. The colors in (b),(d) 
encode the values of the displacement from equilibrium x n 
whereas in (c),(e) the values of the local energies E n defined 
in eq. (18). The figures (b)-(e) show simulation results for 
N = 60 particles and v — 1/2. 


Let us demonstrate this fact with a simple exam¬ 
ple of a Coulomb system consisting of three identical 
charges moving on a ring (Fig. [T^a)). The equilib¬ 
rium configuration of such as system is the equidistant 
one with the charges sitting in the corners of an equi¬ 
lateral triangle with a side of length y. Evidently then 
= y? n = 1,2, 3. If particle 1 is clockwise displaced 
from equilibrium up to the point 1' the resulting local 
excitation energies would be 


energy E n so that it contains: all the parts of the poten¬ 
tial depending solely on u n , none of the parts that are 
independent of u n and half of the terms involving both 
u n and any other U{. This can be written in a compact 
form 

1 N 

E n = 2 E [^(“">“0 + W( Un ,U^)] 

i,i^n 

1 N 

- 2 E mu£\*) + W(u<?\u\ 0) )\+K n , (15) 

i,i^n 


which can be checked to sum up to the total excitation 
energy, i.e. X^=i E n = E. Such a definition solves 
the problem of nonlocality since it distributes all single¬ 
particle potential terms to the corresponding particles, 
yielding for example E n = 0 for u n = Un ^ and K n = 0. 
It can also be checked to solve the problem of negativity 
in many cases, such as the one of the above three-particle 
example. There are however still special cases in which 
negative local energies can be encountered, where the 
negativity originates from the harmonic part of the total 
potential contained in each E n 


y( 2 ) = 

v n 


1 cPV_ 

2 dul 



i,i^n 


(0) V 

d 2 v 

du n dui 



(o) 




which can assume negative values at a given n. To fix 
this, one can make use of the overall stability of the equi¬ 
librium which leads to the existence of a different decom¬ 
position of the harmonic contributions to the potential 
into purely positive parts in the form [28] 

, ( 17 ) 


where Q denotes the square root of the Hessian. We thus 
adjust our definition (15) by subtracting and adding 


instead Vn 2 \ i.e. 


Et 


1 1 

2 y+ + 2 v- 




— e 3 = — -— 

2y' 2 y + 2 y' 


where y+ > y > y_. Obviously E\, E^ > 0, whereas 
E% < 0, demonstrating that following this definition the 
local energy can become negative. Another minor disad¬ 
vantage of this definition of local energy for long range 
interacting particles is that it is highly non-local, mean¬ 
ing that if a single particle is displaced from equilibrium, 
the local energies of other particles even those far apart 
will acquire a finite value. 

A way to fix this last issue is to construct the local 


1 

E n = ^ E f W(Un,Ui) + 

i,i^n 

i,i^n 

+ K n -VW+V™. (18) 

This definition provides a sufficiently local (although 
not entirely, due to the term F n ), positive decom¬ 
position of the excitation energy which still satisfies 
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^2n=i En = E. Its use can be further justified by com¬ 
paring the time evolution of E n with the time evolu¬ 
tion of the displacement of the particle n from equilib¬ 
rium x n = s n — A comparison of Figs. 0b), (d) 

with [7^c) , (e) respectively, leads to the conclusion that 
the time evolution of E n captures nicely all the features 
of the propagating excitation such as the focusing and 
the spreading and only filters out the fast oscillations in 
Figs. 0b) ,(d) caused by the fast continuous conversion 
of potential to kinetic energy and vice versa. 
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